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Automatic, nondestructive methods have been developed for mea- 
suring and analyzing the refractive index distribution of a fiber that 
is immersed in index-matching oil and illuminated transversely to 
its axis in a single-pass interference microscope. The output field of 
the microscope is automatically processed with a video- digitized, 
computer- controlled system, and the profile is determined by the 
solution of an integral equation that can handle arbitrary variations 
in the index distribution. The resulting profiles are reproducible to 
approximately 1 percent and can be determined within a few minutes 
after fiber fabrication. Details of a rapid video scanning procedure 
and of error estimates involved in solving the integral equation are 
presented along with representative profiles. 



I. INTRODUCTION 

One of the most important parameters in determining the usefulness 
of graded-index multimode optical fibers in high-capacity communi- 
cation systems is their refractive index profile. The closer the profile 
conforms to the required optimum distribution, the greater the result- 
ing bandwidth of the fiber. As demand for high-capacity fibers in- 
creases and as fabrication facilities are optimized to produce them 
dependably, it becomes very important to have reliable, accurate, and 
fast methods to measure the fiber's refractive index distribution. These 
methods would greatly aid in tailoring more ideal profiles and would 
form the basis of a fiber evaluation and sorting scheme in which fibers 
are selected for specific applications, depending upon the quality of 
their profiles, almost immediately upon production. 

One of the most sensitive means of measuring the index profile is by 
interference microscopy utilizing cut and polished fiber-slab samples. 1,2 
While this technique allows very precise measurements, 3,4 it requires 
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a high precision, time-consuming procedure of preparing the thin (on 
the order of 50 jum) and extremely flat samples. This can be alleviated 
somewhat by simultaneously processing many samples and utilizing 
an automatic analysis scheme to save time in the evaluation; 5,6 how- 
ever, it would still be preferable to avoid any sample preparation. 

A technique requiring no sample preparation based on interferom- 
etry has been proposed in which a sample of the whole fiber is 
immersed in index-matching oil and illuminated perpendicular to its 
axis. 7 The refractive index profile is subsequently obtained from the 
fringe shifts by mathematical methods which were shown to operate 
for fibers whose profiles are a quadratic function of the radius of the 
core. This method has also been extended 8 to a more general class of 
profile shapes describable by a parameter a which appears in the 
expression describing the index distribution n(r), 



n(r) = /io[l - 2A(r/a) a ] 



1/2 



where r is the distance from the center of the core, a is the core radius, 
and An ■ n — n a . Assuming that a is constant, this treatment showed 
that An and a could be determined to better than ±10 percent. The 
assumption of constant a, however, is a severe restriction not met by 
most currently fabricated fibers and the resulting accuracy of ±10 
percent is also not sufficient for many applications. In addition, these 
methods cannot accommodate arbitrary perturbations of the profile 
which commonly exist in the form of barrier layers at the core-cladding 
interface and an index depression along the axis of the fiber. 

These limitations have been overcome with a method of analysis in 
which profiles with arbitrary shapes can be handled. 9 In that analysis, 
the fiber core is modeled as having many layers, with each layer having 
a constant index of refraction. The profile is then built up in a stepwise 
manner using the index value from each preceding section as the 
starting point for the next. The resulting accuracy is on the order of a 
few percent. 

That analysis, however, is based on having a certain minimum width 
for each layer and breaks down if smaller intervals are utilized. While 
this is not a serious limitation in practice, it may present problems if 
very rapid and narrow variations occur in the profile. We have there- 
fore developed a new analysis scheme not having this potential draw- 
back in which the profile is obtained by the solution of an integral 
equation which can handle small regions. However, this new method 
is still restricted to index profiles of rotational symmetry. We have also 
fully automated the measurement and analysis procedure so that its 
inherent rapid nature can be fully exploited. We call our nondestructive 
measuring procedure the "whole-fiber method." 

This paper consists of two basic sections. In the first, the details of 
the solution of the integral equation and error estimates are treated. 
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In the second, the automated system is described and various profiles 
are presented, discussed, and compared with those obtained by other 
methods. 



II. SOLUTION OF THE INTEGRAL EQUATION AND ERROR ESTIMATE 

In the whole-fiber method, the ray traverses the fiber core at right 
angles to its axis. 8 The phase shift of each ray is an integral over the 
product of the length element of the ray path ds times the relative 
refractive index (relative with respect to the cladding index n c ) . Ignor- 
ing ray bending in the core, we may express the relative phase shift 
according to Fig. 1 as 



w-fjwrt-.j^. 



(1) 



The upper limit of the integral actually need only be p = a but, because 
n(p) = n c for p > a, we may as well use the limit oo. To obtain the 
relation between the fringe shift S(r) and the relative phase shift P, 
we observe that the ratio of fringe shift to fringe spacing D must be 
equal to the ratio of the relative phase shift to 27r, 



S(r) P(r) 



D 



lit 



(2) 




K a H 

Fig. 1 — Outline of fiber core and variables appearing in eq. (1). 
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Equations (1) and (2) thus lead to the following integral equation for 

n(p) - n c \ 



f 



[»(p)-«.]^==4s(r). (3) 

y/p^7 2 2D 



We show in the appendix that the solution of this integral equation 
may be expressed as 



An(r)-n(r)-»--A£ °»r g , (4 ) 



dS(p) dp_ 

dp JJ—p 



This explicit solution of the integral equation (3) shows that the 
refractive index difference at any point r in the fiber core can be 
obtained from a knowledge of the fringe shift function S(r) by differ- 
entiation and integration. Since the fringe shift is known only at certain 
discrete points, numerical techniques [see eq. (14)] for approximating 
the derivative as well as the integral must be used.* The accuracy of 
the resulting refractive index distribution depends on the accuracy of 
the measurement of S(r), on the density of points at which S(r) is 
being measured, and on the sophistication of the methods used for 
numerical evaluation. 

The explicit solution (4) of the integral equation (3) is useful not 
only for obtaining numerical solutions for the refractive index profile 
but also for calculating estimates of the error of An that is introduced 
by the fact that the measured values of the fringe shift S(r) are only 
known to a limited precision. We use the subscript o to indicate the 
index difference A/i that would be obtained if the fringe shift S a and 
its derivative were known precisely. The variance of the refractive 
index difference is thus, according to (4), 

<(An - A/i„) 2 > 

, (a.) 2 r d0 r dD . <*(p> - &o»][*<p'> - & <p-)]> , (5) 

The derivative of S is indicated by the notation S and the symbol 
< > designates an ensemble average. The numerator under the inte- 
gral sign in (5) is the autocorrelation function of the derivative of the 
fringe shift S. For simplicity, we assume that this function is of the 
following form: 



* At and near p = r, an analytical approximation of (4) is used, based on a fourth- 
order polynomial expansion of S(p) . 
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R d = <[£(p) - $ (p)][$( P ') - &(p')]> 

f(A£) 2 for |p-p'|<B 



I 



I R (6) 

p - p \> B. 



The autocorrelation function is assumed to have a constant value, 
equal to the variance of the function S, over a narrow region outside 
of which it vanishes. The interval B assumes the physical meaning of 
a correlation length. We consider B as a very short distance, much 
shorter than any of the transverse dimensions of the fiber. This 
assumption allows us to work out the following approximation for (5). 

«4»-A«.) a >-(^)Wfta(|). (7) 

Approximation (7) is valid only for r » B; however, if we try to 
evaluate (5) with the special function (6) directly at r = we find an 
infinite result. This failure indicates that it is not permissible to use 
(6) for arbitrarily small values of p because the derivative of S, and 
hence its variance, must vanish at p = because of the symmetry of 
the structure. Assuming (A5) 2 to be constant is thus not permissible 
for arbitrarily small values of p. 

The variance (A3) 2 of the derivative of the fringe shift function S is 
not easy to visualize intuitively. For this reason, we must relate it to 
the variance (AS) 2 of the fringe shift function itself. The autocorrela- 
tion function (6) for £ and the autocorrelation function for S 

R(u) = <[S(p) - S (p)][S(p + u) - S„(p + u)]> (8) 

are related by the well-known equation" 

«---!?. (9) 

du 

If we exclude the point | u \ = B, and keep in mind that R ( u ) = R (—u ) 
is required, we conclude that if Rd has the form (6) R must have the 
functional form 



R(u) = 



(AS) 2 B2 B2 U for |m|<£ (1Q) 

for I u I > B. 



Equations (6), (9), and (10) lead to the following relation between the 
variance (A&) 2 of the derivative of the fringe shift function and the 
variance (AS) 2 of the fringe shift function itself: 

2(AS) 2 
~B 



(AS) 2 = ^. (ID 
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„. . ... >/2AAS 
5[An(r)]- ^ 




1/2 


with the definition 


S(An) = [<(An - A«o) 2 >] 


/2 . 



Substituting (11) into the square root of (7) finally yields the following 
expression for the rms deviation of the refractive index difference: 



(12) 



(13) 

It remains to obtain an estimate for the correlation length B. We 
assume that the measured values of the fringe shift S at equidistantly 
spaced measuring points r, are used to define a Lagrange interpolation 
polynomial 

n+l n+l / _ \ 

F= I 11(7-7)* (14) 

»*« 

where we used the abbreviation 

Si = S(n). (15) 

An interpolation polynomial of this type with n = 4 is actually used to 
evaluate the derivative 5 and the integral in (4). We assume that the 
differences of the measured values S, and the true values S (n) are 
mutually uncorrected, 

<(AS,)(AS,)> = (AS,) 2 *,, (16) 

where, by definition, 

AS, = St - S (n). (17) 

We use n+l adjacent points r, to define F according to (14) and 
extend this function outside this range by defining a new function with 
the help of the next set of n + 1 adjacent points, etc. Using 

Rf(u) = <AF(r)AF(r + u)>, (18) 

where AF is defined in analogy with (17), 

&F=F(r)-F (r), (19) 

we have computed the function (18) numerically by using (14) and 
(16). Figures 2 through 5 show the results of these computations for a 
few points r = 1, 1.5, 2, and 2.5 for a fourth-order polynomial approx- 
imation using points n = 1 through r 5 = 5 and extend this function by 
corresponding polynomials outside of this range. The endpoints of one 
range are simultaneously also endpoints of the polynomial expansions 
in adjacent ranges. It is apparent that the autocorrelation function is 
neither stationary nor symmetric because it depends on the choice of 
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Fig. 2 — Correlation function R(r' - r) associated with a fourth-order interpolation 
polynomial if the functional values at the five points used to define the interpolation 
polynomial are mutually uncorrected. The distance between adjacent interpolation 
points is used as the unit of length. In this figure, the autocorrelation function is 
evaluated with r = 1 being the end point of adjacent interpolation intervals. Adjacent 
intervals (of 5 points each) are approximated by similar interpolation polynomials. Note 
that, for this choice of r, the function assumes a symmetrical appearance. 



r in (18). On the other hand, we may use this approach to gain an 
order-of-magnitude estimate of the correlation length B. Figure 6 
shows an average of the autocorrelation functions averaged over r. 
This function is very nearly of the form (sin x)/x and, in its central 
portion, is a reasonable order-of-magnitude approximation of the cor- 
relation function (10). It is clear from this averaged autocorrelation 
function that the correlation length B is on the order of the length of 
the interval used for defining the individual interpolation function, 

B = nAr, (20) 

where Ar indicates the distance between adjacent sampling points 

Ar = r 1+1 - n. (21) 

Even though our error analysis can claim to give no more than an 
order-of-magnitude estimate, (12) and (20) show that the error due to 
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Fig. 3— Similar to Fig. 2 with r = 1.5. The autocorrelation function is now not 
symmetrical around the point r' - r = 0. 

the limited accuracy of the measurements of the fringe shift S(r) 
increases as r decreases and also that the error of the refractive index 
determination increases as the length of the interval Ar between 
adjacent sampling points decreases. This latter error can be decreased 
by using an interpolation polynomial of higher order. All these trends 
have been clearly observed in our experimental results. 

To gain insight into the numerical values of the rms deviation 8(An), 
let us make the following assumptions. The distance between adjacent 
scan lines of the vidicon corresponds to 0.2 tun of the vertical direction 
along the fiber core. Because of an averaging procedure using 10 
adjacent scan lines for finding the center of each fringe, we may assume 
AS = 0.2/10 = 0.02 /un. The distance between adjacent sampling points 
in horizontal direction is typically Ar = 1 /xm. We set n = 4 because of 
our use of a fourth-order interpolation polynomial. The fringe spacing 
is typically D = 6 iim. If we use A = 1 iim and use a value halfway 
between core center and core boundary, r = 15 /im, we find, from (12), 



5(An) = 3 X 10" 



(22) 



This amount of rms deviation is in excellent agreement with the 
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Fig. 4— Similar to Figs. 2 and 3 with r = 2. 



observed values. However, it is important to remind the reader that 
this is the error caused by the random uncertainty in the fringe position 
that is inherent in the measuring process. It is not the absolute error 
relative to the (unknown) precise value of An. Our derivation of the 
error estimate (12) is limited to the random component of the uncer- 
tainty of the fringe shift S and does not include systematic errors 
arising from distortions in the imaging process of the microscope or 
the vidicon and furthermore does not include the systematic error 
arising from the process of using a finite number of discrete points for 
evaluating the derivative 5 and the integral in (4). Only the additional 
error caused by the uncertainty about the actual values of S(n) at the 
sampling points is responsible for the rms deviation of An expressed 
by (12). Systematic errors existing in the system have been investi- 
gated, however, and were found to be at least one order of magnitude 
less than this value. 

III. AUTOMATIC PROCESSING AND ANALYSIS SYSTEM 

The experimental arrangement to automatically measure and ana- 
lyze the profiles of whole fibers is shown in Fig. 7. The use of this 
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Fig. 5 — Similar to Figs. 2 and 3 with r = 2.5. 



system has been briefly described elsewhere. 9 Here, we discuss the full 
details of its operation. 

The basis of the system is a Leitz dual-beam, single-pass transmis- 
sion interference microscope which is generally used to perform precise 
refractive index profiling of optical fibers by examining polished slab 
samples. 1 " 6 For the whole-fiber measurements, a short length of fiber 
(~1 cm) is inserted into index-matching oil in the sample arm of the 
microscope and a similar thickness of oil is placed in the reference 
beam. An excellent match to the cladding could be achieved at an 
observation wavelength \ = 0.9 fim with matching oil of index i\d = 
1.457 ± 0.0005 with no temperature control. 

The output field of the microscope is detected with an infrared- 
enhanced, silicon-target vidicon, and the video signal is sent to a video 
digitizer that has the capability of addressing discrete picture elements 
in the television frame and digitally encoding the intensity at each 
element. The digitizer resolves 480 picture elements along the Y axis 
and 512 elements on the X axis. The X- and Y-position data inputs are 
provided by the 16-bit duplex input/output (I/O) interfaces of a 
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Fig. 6— Averaged autocorrelation function. The average is taken over r in the interval 

1 < r == 5. 
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Fig. 7— Experimental arrangement to perform automatic index profiling of whole 
fibers. 
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Hewlett Packard 9825A computer with the least significant bit of the 
Y address serving as the field selector. 

The digitizer normally operates by deriving horizontal (line) and 
vertical (field) timing from the video input signal and uses these timing 
signals to keep track of the position of the video signal relative to the 
TV raster at any instant. This location is compared with the selected 
X, Y sample point and, when coincidence is achieved, the video at that 
point is digitally encoded to 8 bits or 256 gray levels and sent to the 
computer which then issues new X, Y coordinates and the procedure 
is repeated. We have found, however, that the transition of the com- 
puter between issuing new data (writing-mode) and receiving the 
digitized information (reading- mode) is sufficiently long (~3 ms) that 
it does not allow for the digitization of all Y values for a given X 
coordinate (one column) within one field, 1/60 s. To do this requires 
an execution time of 50 fis or less. In fact, only one point per field can 
be processed, and digitizing a single column typically takes about 3 s. 
This is not normally a serious limitation since only on the order of 50 
columns (requiring 2Vfe minutes) are needed to measure the profile with 
good resolution. However, if it is desired to obtain maximum sensitivity 
by averaging over many scans to reduce video noise levels, this limi- 
tation becomes prohibitive. 

The computer, on the other hand, is capable of reading 400,000 
words of data per second and if it could be freed of the need to supply 
consecutive Y addresses, digitization of a column could be accom- 
plished in 1/60 s (for each field), saving a factor of almost 200 in time. 
This has been achieved with the addition of the interface board 
shown schematically in Fig. 8, which essentially contains an 8-bit 
programmable counter and support logic gates. At the start of each 
columnar digitization, the computer loads the counter with the logic 
value of the number of points desired and issues the initial X value. 
These readings are then taken, and a stop pulse is issued. An enable 
pulse whose length is related to the number of data points is also 
generated by the circuitry to supply handshaking (timing) for the 
system and to allow the region being processed to be observed on a 
display monitor. 

A photograph of the display showing a typical graded-index fiber 
observed at a wavelength of A = 0.9 /un is seen in Fig. 9. Location of 
the point being encoded is indicated by a dot cursor and the gray scale 
value for all Y elements along a sample line at the selected X position 
is displayed as a waveform at the side of the monitor screen. The two 
vertical lines bounding the waveform display indicate the black and 
white reference levels which are encoded as 255 and 0, respectively. 

The sequence of events for determining the fringe shift as a function 
of radius, S(r), from which the profile is computed, shall be described 
with reference to Fig. 9. It is first necessary to choose a desired fringe. 
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Fig. 9 — Display of video monitor showing transversely illuminated whole fiber at a 
wavelength of 0.9 jun. 



This generally is not critical since any fringe should give the same 
index distribution assuming that the profile is constant along the fiber 
and the sample is relatively clean and free of surface imperfections. 
Once selected, the fringe is located by a procedure of fringe counting. 
For example, starting from point A and proceeding along a line AB, let 
the desired fringe number be 3. Thus all elements along line AB 
defined by a fixed X value are digitized. The computer then counts the 
fringes to find the third one. In this procedure, a fringe is defined by 
following the digitized data and noting transitions in gray-scale values 
by more than some fixed number. Twenty data points in the vicinity 
of the fringe maxima are chosen, and a parabolic fit to them is formed. 
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The maximum of the curve thus determined is taken as the center of 
the fringe. The computer also determines at this time the uniform 
spacing, D, of the fringes in the cladding by locating the center of each 
fringe and averaging over the number of fringes. 

The X position value is then changed and the same determination 
is made for the right-hand side of the cladding along line CD. A 
straight line is fitted by the computer to pass through the centers of 
the third fringe at these right- and left-hand sides. This line serves as 
the reference level from which fringe displacements are made and also 
compensates for possible misalignment of the fringes to the scanning 
line axis. 

Next, the X position is automatically set to point E which is about 
midway between fringe 3 and 4, and scanning proceeds from EtoB.X 
is then incremented for each such columnar encoding and the Y 
address is automatically corrected to track the fringe by picking 
midway initial X values as was done for point E. 

To enable the calculation of the radius, the XY coordinates of the 
fiber's axis are located as part of an initial set-up procedure. This is 
readily achievable due to the index depression which exists along the 
axis 10 and causes a perceptible dip in the fringes. The initial set-up 
procedure also includes entering into the computer the coordinates of 
points A, B, and C, the wavelength of observation, the desired X 
increment, and the desired fringe number. The computer then takes 
over and at the end of the encoding process calculates An as a function 
of radius for each half of the core, averages both halves, and plots the 
resulting index profile along with coordinates and labeling on an XY 
plotter. 

A best-fit power-law (a) curve to the index profile may be determined 
after the coordinates of the limits over which the fit is to be made are 
entered. In addition to the a parameter, this fitting procedure also 
determines the core radius, the shift of the center of the core, the 
maximum index difference between core and cladding, and a fitting 
error to provide a measure of how well the a curve approximates the 
measured profile. These values are then printed out and the best-fit 
a curve is drawn on the plotter. 

Four examples of measured profiles for various Ge0 2 -doped fibers 
demonstrating the general applicability of this method are shown in 
Figs. 10 to 13. 

Figures 10 and 11 show typical and reasonably smooth profiles (solid 
line) which are fit with relatively good a curves (broken lines), the 
fitting errors being on the order of 1 percent. Perturbations in the 
profile, mainly in the form of the index depression along the axis and 
variations in deposition layers close to the center, are clearly resolved. 
It is important to emphasize that this method assumes circular sym- 
metry for the fiber core and thus variations in core geometry, if 
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Fig. 10— Profile of typical Ge0 2 -doped fiber (solid curve) and best-fit a curve (broken 
curve). 

existing, would not be detected. The variations in the plotted profile 
also represent averages over both sides of the core. 

Figure 12 shows the profile of a fiber with strong index variations in 
the vicinity of the axis, and Fig. 13 shows a fiber with a very large 
boron barrier layer at the core-cladding interface. These unique fea- 
tures are clearly detectable and well resolved, thus demonstrating the 
applicability of this technique to even these extreme cases. 

We have compared profiles obtained using the circular index analysis 
technique of Ref. 8, in which the profile is built up piecewise, with the 
integral equation approach used here. The agreement was excellent 
with the profiles being essentially identical. 

The profiles have also been compared with those obtained by 
preparing polished slabs of the same fiber samples. The comparison 
results, presented in detail in Ref. 8, show very good agreement with 
the average difference of the maximum An values and the a values 
being about 4 percent. 
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The experimental accuracy of the measured profiles and the effects 
of index mismatching were investigated by measurements made on 
unclad fiber (uniform- index) samples under various matched and un- 
matched conditions. Variations of the index were consistently found 
to be within about 2 parts in 10 4 , in excellent agreement with the 
theoretical analysis of the previous section. For unclad fibers with 
intentional index mismatch, we obtained step index profiles with the 
height of the step corresponding to the degree of mismatch. When a 
good match is visually observed on the monitor, as evidenced by a 
straight fringe parallel to the scanning line axis, the error in the profile 
due to mismatch was found to be on the same order as that due to the 
index variations (a few parts in 10 4 ); thus, no attempt at more precise 
matching was made. 

In conclusion, the techniques described to measure the refractive 
index profiles of optical fibers rapidly and automatically should prove 
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Fig. 11— Profile of Ge0 2 -doped fiber (solid curve) and best-fit a curve (broken curve) 
showing somewhat greater variations than fiber of Fig. 10. 
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Fig. 12 — Profile of GeC>2-doped fiber (solid curve) and best-fit a curve (broken curve) 
showing strong central index perturbations. 



valuable not only in reducing the time required for these measurements 
but also in forming the basis of an evaluation scheme in which fibers 
can be rated and selected for specific applications, depending upon the 
quality of their profiles, almost immediately upon production. 



APPENDIX 

We sketch the solution of the integral equation (3), which we write 
in the form 



xf(x) dx 
yy/(x 2 -/) 



= g(y)> 



(23) 



In our application, f(x) is an even function of x; this allows us to 
express it in the form 



fix) 



= J F(u 

Jo 



)cos(ux) du. 



(24) 



Substitution of (24) into (23) yields after integration 

900 THE BELL SYSTEM TECHNICAL JOURNAL, APRIL 1979 



-I F(u)Ji(uy)du = -g{y), 



(25) 



where J\ is the Bessel function of order one. This integral relation is 
the Hankel transform whose inversion is given by the formula 



F(u) = --u giy)Ji(uy) dy. 



(26) 



Substitution in (24) yields 



f(x) = - - | g(y)\ | ucos(ux) Ji(uy) du 



\dy. 



(27) 



19A118A 
in = 0.0088 
r = 24.47/Xm 
a = 2.15 
FIT = 2.36% 




Fig. 13— Profile of Ge02-doped fiber (solid curve) and best-fit a curve (broken curve) 
showing resolution of large boron barrier layer. 
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This expression can also be written in the form: 

fix) - 1 j g(y)lj- [ cos(ux)J (uy) du\ dy. (28) 

We perform a partial integration 
fix) 

= --\g(0) cos(ux)du+ I ~ I cos(ux)J (uy) dudy\, (29) 



7T 



J —I cos(ux)J (uy) dudy>, 



where we used the fact that the fringe shift vanishes outside the fiber 
core so that g(<*>) = 0. The integral over the cosine function vanishes 
for x ¥* 0. The integral over the product of the cosine and Bessel 
functions can be found in tables so that we finally obtain 






dg dy 
dy Jy*—?' 



The integral transform pair (23) and (30) corresponds to the transform 
pair (3) and (4). 
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